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Abstract 

In the reduced basis method, the evaluation of the a posteriori estimator 
can become very sensitive to round-off errors. In this note, the origin of the 
loss of accuracy is revealed, and a solution to this problem is proposed and 
illustrated on a simple example. 

1 Introduction 

Consider the following discrete variational form, depending on a parame- 
ter \i e V : find in a Hilbert space V such that Vt> G V, : a p {u^v) = 
b(v), where (X^ is a bilinear form and b a linear form. In a many queries 
context, a quantity of interest of the solution q{u^) has to be computed for 
many values of \i E V . In this note, the following assumptions are made for 
simplicity, but the conclusions are general : (i) the variational formulation 
is coercive, (ii) has a so-called affine dependence on the parameter fi so 
that a M = a + ai(/i)a 1; (iii) the quantity of interest is the solution itself : 
g(it M ) = v 

Reduced basis (RB) strategies consist in replacing E^ by an easily com- 
putable surrogate E^, that is made precise below (see [2]). We denote 
the solution of E^, N the size of the matrix involved in the resolution of 
E^, and N the size of the matrix involved in the resolution of E^. The RB 
method consists in two steps : (i) An offline stage, where a basis, whose 
vectors are solutions of E^ for well-chosen values of the parameter /i, is 
constructed using, e.g., a greedy algorithm on the parameter. During this 
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stage, N problems of size N are solved, and some quantities related to the 
solutions are stored, (ii) An online stage, where the precomputed quantities 
are used to solve for many values of fi. In this stage, an a posteriori 
error estimator £(/x) is also computed to check the quality of the approxi- 
mation. This is called certification. The a posteriori error estimator verifies 
\\v>n — u n\\v ^ £(/-0 := /^II^m^vIIvj where /3 M is the coercivity constant of 
(or a lower bound of it) and is the unique afline application from 
V to V such that \/(u,v) G V 2 , (Gyu,t>)y = a^(u,v) — b(v). In this note, 
we consider different ways to compute the same quantity £(//). We distin- 
guish between formulae to compute £(//) by adding an index to £(//). Thus, 
£i(n) := /3 A r 1 ||G r M 'U M ||v is the first formula for the estimator, directly given 
by the definition. Since £i(fj) requires the computation of a size iV scalar 
product, this formula is not compatible with the constraint that the compu- 
tations in the online stage should be of complexity independent of N. 

Suppose that a reduced basis of size N has been computed in the offline 
stage, namely a family Ui) i=1 where each of the iV-dimensional vectors 
Ui is the solutions of E^.. The reduced problem is a Galerkin procedure on 
the Ui basis : find w M G Span{ui, u^} such that Vj G {1,...,N}, E^ : 
a^u^Uj) = b(uj). Writing ii^ = Yld=ili{^) u i-> the reduced linear system 
is AJjf, = b, where (A^ij = a^u^Uj), (b)j = b(uj) and (U^i = 7i(/i). 
Using the afline parameter dependance, the matrix (A^)^ = a (ui,Uj) + 
ai(/i)ai(ui,Uj) is built and solved in complexity independent of N, provided 
that the quantities a (ui,Uj) and ai(ui,Uj) have been precomputed in the 
offline stage. 

Consider the Riesz isomorphism J from V' to V such that V/ G V, Vm G V, 
(Jl,u) v = l(u). The operator G M inherits the affine dependance of a M in /i 
since, Vw G V, 

Gyu = —Jb(-) + Ja (u, •) + Ja^u, •) =: G o + G w + ot\{n)GiU, (1) 

where &(•) and a k (u, •), G {0, 1}, are elements of V. The a posteriori error 
estimator is then written in the following compact form : 

£ 2 (/i) := (5 2 + 28% + x^Sx,) * , (2) 

where 5 = ||G o||v, = (Goo, G k Ui) v , x M = a fc (/i)7i(/i), S/,j = (G k u u GiUj) v 
(with / and J re-indexing respectively (fc, i) and (/, j), ^ A;, Z ^ 1, 1 ^ i,j ^ 
iV) and a = 1. Provided that 5 G R, s G M 2i \ and 5 G ]R 2 * x2 ^ (which 
are independent of //) have been precomputed in the offline stage, ^(/-O is 
computed in complexity independent of N. This is what is typically used in 
RB implementations. 
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2 Round-off errors and certification 



Canuto, Tonn and Urban [T] identified that the evaluation of ^(aO suffers 
in practice from a loss of accuracy, which they attributed to the square root 
in [2j Herein, we show more precisely that this loss of accuracy comes from 
round-off errors. Indeed, when substracting two real numbers within floating 
point arithmetics, the number of lost significant digits equals the number 
of common decimals between the two reals. For simplicity, we neglect the 
round-off errors introduced when solving and E^, so that the vectors of 
the reduced basis U{ and the reduced solutions are considered free of round- 
off errors. Therefore, we only consider round-off errors in the evaluation of 
and ^(aO due to the summations. We define the machine precision e 
by the maximal floating point representation relative error of real numbers : 
' /<(«)- 



^ e. Under these hypotheses, the smallest possible values that can 

3e practically computed for Si(fi) and ^(/-O using floating point arithmetics 
with machine precison e is bounded below by respectively 4-e and -i-y/s. This 
is supported numerically (see section [4]). 

This observation is of paramount importance since the certification of a 
RB procedure cannot be better than these values. In a successful RB pro- 
cedure, the value of the estimator gets smaller as the size N of the reduced 
basis increases. Enriching the basis with a new vector improves the quality 
of the approximation introduced by the method. As a result, there exists N 
such that, ViV > A> , v> G V, £(//) < -f^. If N > N Q , £ 2 (n) is no longer 
suitable for computing the a posteriori error estimator. 

We notice that increasing the machine precision from e to e 2 enables the 
accuracy of ^(aO to reach the one of Thus, the use of quadruple 

precision is a first solution, checked numerically in section |4j This is however 
not pratical since current computer architectures are optimized for double 
precision. Another solution is to develop an alternative algorithm for the 
evaluation of the estimator that still achieves the machine precision. This is 
the purpose of the next section. 



3 The new procedure for a posteriori error 
evaluation 

Consider that a reduced basis of size N has been constructed. Let us 
denote d = 1 + 3N + 2N 2 . For a given \i and G Span{wi, u^}, we define 
X(fi) G M d as the vector with components (1, x^j, x^jX^j), with 1 ^ / ^ 
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J ^ 2N. Using the symmetry of the matrix S, we can write the right-hand 
side of 2 as a linear form in X(fi) : J2p=i QpXpifi), where q p is independent 
of n ana X p (/i) is the p-th component of X(fi). 

During the offline stage, we take d values, possibly random, /i r , r £ 
{1, d} of the parameter \i. Then, we compute the vectors X(fj, r ) and, using 
the accurate formula E\ for the estimator, the quantities V r := (/3 Mr £i(/z r )) 2 . 
Finally, we define T £ M> dxd as the matrix whose columns are formed by the 
vectors X(fi r ) and we assume that T is invertible, which was the case in our 
simulations. 

In the online stage, suppose that we want to evaluate the estimator for the 
parameter value fi. We compute the vector X(/a) and solve the linear system 
T\(ji) = X(ji), for AO) £ R d . We then have X(ji) = Y? r =i Kip) x M and 

d d d d 

p=l p,r=l r=l r=l 

; (3) 

This yields the new formula for computing the estimator, 



:= ^ ( X>MV r ] . (4) 

Quite importantly, we notice that the additional cost is such that the quantity 
B^ip) is still computed in complexity independent of N. 

The quantity (^^(/x)) 2 is a sum of terms, whose first one is fixed and 
equals 5 2 . On the contrary, ((3^S 3 (n)) 2 is a sum whose terms V r are updated 
each time a vector is added to the reduced basis. Since V r is computed for /x r £ 
V, (f3fj,£(fi)) 2 and V r , (or, at least, max^gp (/3 M £^(/i)) 2 and max^^,! V r , which 
are the important quantities in a greedy selection and online certification) 
are of the same orders of magnitude. 



4 Numerical illustration 

Consider the equation —u" + fiu = 1 on ]0, 1[ with u(0) = u(l) = 0, where 
li ^ 1 is the parameter. The analytic solution is u(x) = — - (cosh (y/jlx) — l) + 

cosh sm h (y//ix) . The Lax-Milgram theory is valid, the coercivity constant 

is 1 in the if 1 -norm. The estimator is given by £((jl) = ||G>^||.e- 1 (]o 1 in- La- 
grange Pi finite elements are used, with uniform mesh cells of length 0.005. 
The RB method is carried-out until a reduced basis of size 6 is constructed. 
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Figure 1 - Left : £%, £3 algorithms for the estimator and error curves with 
respect to the parameter //; Right : £1, £2 (double and quadruple precision) 
curves. 



On the left part of Figure [TJ £ 3 yields the same curve as the accurate 
but expensive £\ algorithm. Notice that the values of the a posteriori error 
estimators are very close to the values of the error. This means that the 
efficiency of the estimator is very close to 1. On the right part of Figure 
[1} we see that £2 yields a flat curve for the estimator, meaning that all 
information relative to the error is lost. As expected, the use of quadruple 
precision enables £2 to recover the accuracy levels of £\. In this example, 
j rj 3 x 10~ 17 and ^ ^3x 1CT 9 , which sould be respectively compared to 
the numerical values of £ 2 in quadruple precision and £\ (1CT 17 and 10~ 15 ) 
and £2 (1IT 8 ). 



5 Conclusion 

To sum up, we have developed a procedure where the accuracy of the 
online evaluation is limited by the accuracy of the evaluation of quantities 
precomputed during the offline stage, where heavy but accurate algorithms 
are allowed. In the online stage, instead of a linear combination of d terms, we 
have to solve a linear system of size d, before doing a linear combination of the 
same size. We have increased the accuracy of the estimator, with a procedure 
of complexity independent of the size N of the initial problem. When the size 
of the reduced basis increases, we observe that the condition number of the 
matrix T increases as well. Finally, we notice that oversampling strategies 
consisting in defining a least squares problem to compute A(/x) such that T 
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is rectangular with more than d columns improve the quality of our results 
when the RB is close to convergence. Experiments on a more complicated 
problem (external acoustics solved by integral equations where the criterion 
is on the far field approximation of the diffracted acoustic potential) lead to 
similar conclusions, which will be discussed in more detail elsewhere. 
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